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Abstract 

The corrections to finite-size scaling in the critical two-point correlation 
function G(r) of 2D Ising model on a square lattice have been studied numer- 
ically by means of exact transfer-matrix algorithms. The systems of square 
geometry with periodic boundaries oriented either along (10) or along (11) di- 
rection have been considered, including up to 800 spins. The calculation of 
G{r) at a distance r equal to the half of the system size L shows the existence 
of an amplitude correction cx L^^. A nontrivial correction cx L""'^^ of a very 
small magnitude also has been detected in agreement with predictions of our 
recently developed GFD (grouping of Feynman diagrams) theory. A refined 
analysis of the recent MC data for 3D Ising, ip"^, and XY lattice models has 
been performed. It includes an analysis of the partition function zeros of 3D 
Ising model, an estimation of the correction-to-scaling exponent lu from the 
Binder cumulant data near criticality, as well as a study of the effective critical 
exponent rj and the effective amplitudes in the asymptotic expansion of suscep- 
tibility at the critical point. In all cases a refined analysis is consistent with our 
(GFD) asymptotic values of the critical exponents (i/ = 2/3, u; = 1/2, rj = 1/8 
for 3D Ising model and w = 5/9 for 3D XY model), while the actually accepted 
"conventional" exponents are, in fact, effective exponents which are valid for 
approximation of the finite-size scaling behavior of not too large systems. 

Keywords: Transfer matrix, Ising model, XY model, model, critical expo- 
nents, finite-size scaling, Monte Carlo simulation 
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1 Introduction 



Since the exact solution of two-dimensional Ising model has been found by On- 
sager [|l|, a study of various phase transition models is of permanent interest. Nowa- 
days, phase transitions and critical phenomena is one of the most widely investigated 
fields of physics P, ^. Remarkable progress has been reached in exact solution of 
two-dimensional models Q. Recently, we have proposed |^] a novel method based 
on grouping of Feynman diagrams (GFD) in (p^ model. Our GFD theory allows to 
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analyze the asymptotic solution for the two-point correlation function at and near 
criticality, not cutting the perturbation series. As a result the possible values of exact 
critical exponents have been obtained |^] for the Ginzburg-Landau (93^) model with 
0{n) symmetry, where n = 1,2,3,... is the dimensionality of the order parameter. 
Our predictions completely agree with the known exact and rigorous results in two 
dimensions Q, and are equally valid also in three dimensions. In we have com- 
pared our results to some Monte Carlo (MC) simulations and experiments 0, ^. 
It has been shown |^] that the actually discussed MC data for 3D Ising |^] and 
XY models are fully consistent with our theoretical predictions, but not with 
those of the perturbative renormalization group (RG) theory [10, 11, 12|. From the 
theoretical and mathematical point of view, the invalidity of the conventional RG 
expansions has been demonstrated in ||5|. The current paper, dealing with numerical 
transfer-matrix analysis of the two-point correlation function in 2D Ising model, as 
well as with the analysis of MC data for the three-dimensional Ising, Xip'^ and XY 
models, presents a more general confirmation that the correct values of critical ex- 
ponents are those predicted by the GFD theory. Our estimations are based on the 
finite-size scaling theory, which by itself is an attractive field of investigations |13] 
and has increasing importance in modern physics 



2 Critical exponents predicted by GFD theory 

Our theory provides possible values of exact critical exponents 7 and u for the (^^ 
model whith 0(n) symmetry (n-component vector model) given by the Hamiltonian 

H/T = J [roip\x) + c(V(^(x))2 + m99^(x)] dx , (1) 

where tq is the only parameter depending on temperature T, and the dependence is 
linear. At the spatial dimensionality d = 2, 3 and n = 1, 2, 3, ... the critical exponents 
are 

d + 2j + 4m 



d{l + m + j)- 2j 

2(1 +m) +i 
d{l +m + j)- 2j 



(3) 



where m > 1 and j > —m are integers. At n = 1 we have m = 3 and j = to fit the 
known exact results for the two-dimensional Ising model. As proposed in Ref. [^, 
in the case of n = 2 we have m = 3 and j = 1, which yields in three dimensions 
u = 9/13 and 7 = 17/13. 

In the present analysis the correction-to-scaling exponent 9 for the susceptibility 
is also relevant. The susceptibility is related to the correlation function in the Fourier 
representation G(k), i. e., x oc G(0) [|ll|. In the thermodynamic limit, this relation 
makes sense at T > Tc, where Tc is the critical temperature. According to our 
theory, t'^G{0) can be expanded in a Taylor series of t^"-'^ at t ^ 0. In this case the 
reduced temperature t is defined as t = ro(T) — ro(Tc) ccT — Tc- Formally, ^27-^!^ 
appears as second expansion parameter in the derivations in Ref. Q, but, according 
to the final result represented by Eqs. (^) and (^, (27 — diy)/{2i/ — 7) is a natural 
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number. Some of the expansion coefficients can be zero, so that in general we have 



e = l{2u-^), (4) 

where £ may have integer values 1, 2, 3, etc. One can expect that i = A holds at 
n = 1 (which yields 9 = 1 at d = 2 and 6 = 1/3 at d = 3) and the only nonva- 
nishing corrections are those of the order , t^^ , t^^, since the known corrections 
to scaling for physical quantities, such as magnetization or correlation length, are 
analytical in the case of the two-dimensional Ising model. Here we suppose that 
the confluent corrections become analytical, i. e. 6 takes the value 1, at d = 2. 
Besides, similar corrections to scaling are expected for susceptibility x mag- 
netization M since both these quantities are related to G(0), i. e., x G{0) and 

= limx^oo{v{0)(p{x)) = limy^oo ^(0)/^ hold where V = L'^ is the volume and 
L is the linear size of the system. The above limit is meaningful at L ^ oo, but 
G{0)/V may be considered as a definition of for finite systems too. The latter 
means that corrections to finite-size scaling for x and M are similar at T = Tc. Ac- 
cording to the scaling hypothesis and finite-size scaling theory, the same is true for 
the discussed here corrections att^O, where in both cases (x and M) the definition 
t =1 ro(r) — rQ{Tc) I is valid. Thus, the expected expansion of the susceptibility x 
looks like x = t^"' (^0 + ait^ + ^2*^^ H ) • 

Our hypothesis is that j = j{n) and i = i{n) monotoneously increase with n to 
fit the known exponents for the spherical model at n ^ oo. The analysis of the MC 
and experimental results here and in enables us to propose that j{n) = n — 1, 
i{n) = n + 3, and m = 3 hold at least at n = 1,2. These relations, probably, are 
true also at n > 3. This general hypothesis is consistent with the idea that the 
critical exponents 7, z^, and 9 can be represented by some analytical functions of n 
which are valid for all natural positive n and yield r] = 2 — ^ /v (x 1/n rather than 
r] oc 1/n* with s = 2,3, . . . {s must be a natural number to avoid a contradiction, 
i. e., irrational values of j{n) at natural n) at n — > 00. At these conditions, j{n) 
and t{n) are linear functions of n (with integer coefficients) such that i{n)/j{n) 1 
at n 00, and m is constant. Besides, j(l) = 0, m(l) = 3, and ^(1) = 4 hold 
to coincide with the known results at n = 1. Then, our specific choice is the best 
one among few possibilities providing more or less reasonable agreement with the 
actually discussed numerical an experimental results. 

We allow that different i values correspond to the leading correction-to-scaling 
exponent for different quantities related to G(k). The expansion of G(k) by itself 
contains a nonvanishing term of order t'^'^^'^ = t'^'^ (in the form G(k) ~ 
t-T [gikt-") + f^" giikt'")] whith 5-1(0) = 0, since £>1 holds in the case of suscep- 
tibility) to compensate the corresponding correction term (produced by c {Vip)"^) in 
the equation for 1/G(k) (cf. Q). The latter means, e. g., that the correlation length 



^ estimated from an approximate ansatz like G(k) oc 1/ [k^ -|- (1/^)'^] used in |14, 15] 
also contains a correction proportional to t'^'^ . Since "qv has a rather small value, 
the presence of such a correction (and, presumably, also the higher order correc- 
tions t'^^'^ , t^^'^ , etc.) makes the above ansatz unsuitable for an accurate numerical 
correction-to-scaling analysis. 

The correction t'^'^ is related to the correction L~'^ in the finite-size scaling expres- 
sions at criticality. The existence of such a correction in the asymptotic expansion 
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of the critical real-space Green's (correlation) function is confirmed by our results 
for the 2D Ising model discussed in Sec. [4.3[ 

Our consideration can be generalized easily to the case where the Hamiltonian 
parameter rg is a nonlinear analytical function of T. Nothing is changed in the 
above expansions if the reduced temperature t, as before, is defined hy t = ro(T) — 
ro(Tc). However, analytical corrections to scaling appear (and also corrections like 
(T — 7"^)™+"-^ with integer m and n) if t is reexpanded in terms of T — Tc at T > Tc- 
The solution at the critical point remains unchanged, since the phase transition 
occurs at the same (critical) value of tq. 

3 Exact transfer matrix algorithms for calculation of 
the correlation function in 2D Ising model 

3.1 Adoption of standard methods 

The transfer matrix method, applied to analytical calculations on two-dimensional 
lattices, is well known ^. However, no analytical methods exist for an exact cal- 
culation of the correlation function in 2D Ising model. This can be done numerically 
by adopting the conventional transfer matrix method and modifying it to reach the 
maximal result (calculation of as far as possible larger system) with minimal number 
of arithmetic operations, as discussed further on. 

We consider the two-dimensional Ising model where spins are located either on 
the lattice of dimensions N x L, illustrated in Fig. |^a, or on the lattice of dimensions 
\/2N X V2L, shown in Fig. ||b. The periodic boundaries are indicated by dashed 
lines. In case (a) we have L rows, and in case (b) - 2L rows, each containing spins. 
Fig. m shows an illustrative example with = 4 and L = 3. Let us cr{k) be the spin 
variable (±1) in the fc-th node of the first row. Here nodes are numbered sequently 
from left to right, and rows - from bottom to top. Our method can be used to 
calculate the correlation (Greens) function between any two spins on the lattice. As 
an example we consider the Greens function G{r) in (10) crystallographyc direction, 
indicated in Fig. ^ by arrows, i. e., 

G{r) = {a{k)(j{k + r)) : case (a) (5) 

G{r) = {a{k)a'{k + A{r))) : case (b) . (6) 

Here a' refers to the (1 + r)-th row. It has a shift in the argument k by A(r) = r/2 
for even and A(r) = (r — l)/2 for odd r. It is supposed that a{k + A^) = cr{k) holds 
according to the periodic boundary conditions. 

For convenience, first we consider an application of the transfer matrix method 
to calculation of the partition function 

Z=J2eKp[(3Y,aiaA , (7) 

Wk] \ (ij) J 

where the summation runs over all the possible spin configurations {(t^}, and the 
argument of the exponent represents the Hamiltonian of the system including sum- 
mation over all the neighbouring spin pairs {i,j) of the given configuration {dk}', 
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Figure 1: Illustrative examples of the lattices with dimensions NxL{si) and V2N x ^/2i (b) 
with periodic boundary conditions along the dashed lines. The correlation function has been 
calculated in the (10) crystallographic direction, as indicated by the arrows. 



parameter P is the coupling constant. Let us consider lattice (a) in Fig. ||, but 
containing n rows without periodic boundaries along the vertical axis and without 
interaction between spins in the upper row. We define the 2'^-component vector 
r„ such that the i-th component of this vector represents the contribution to the 
partition function provided by the i-th spin configuration of the upper row. Then 
we have an obvious recurrence relation 

r„+i = rr„, (8) 

where T is the transfer matrix with the elements 

Tij = exp [^(^ + % + /3E [^(^)].) ■ (9) 

\ fc k ) 

Here [(7(A;)]^ is the spin variable in the /c~th position in a row provided that the 
whole set of spin variables of this row forms the z-th configuration. The first and 
the second sum in (|9|) represent the Boltzmann weights for the spin interaction in 
the ri-th row, and between the n-th and (n + l)-th rows, respectively. In this case, 
ri is the vector with components (ri)j = 1- If we set 

(ri), = , (10) 

then the components of the resulting vector r„ = r„ give us the partial contributions 
to the partition function corresponding to a fixed, i. e. i-th, configuration of the first 
row. The periodic boundary conditions along the vertical axis means that the (L + 
l)-th row must be identical to the first one, i. e., we have to take the i-th component 
of the vector r]^i|_^ and make the summation over i to get the partition function Z 
of the originally defined lattice in Fig. |l|a. Note that the missing Boltzmann weights 
for the interaction between spins in the {L + l)-th row are already included in the 



first row. By virtue of (^) and ([10|), we arrive to the well known expression [Q, 16 1 



^ = E(r?li), = Trace(r^)=EAf , (11) 
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where Aj are the eigenvalues of the transfer matrix T. An analogous expression for 
the lattice in Fig. [11) reads 



Si) 
2L+1 



Trace ([r2ri] 



(12) 



where the vectors r„ obey the reccurence relation 



(13) 



similar to (^), but with different transfer matrices Ti and T2 for odd and even row 
numbers n, respectively. They include the Boltzmann weights for the interaction 
between two subsequent (odd-even or even-odd, respectively) rows, i. e.. 



(Ti,2),, = exp IpYl i<^ik)]^ {Hk)]. + [a{k ± 1)]^.} 



(14) 



The actual scheme can be easily adopted to calculate the correlation functions (|5|) 
and (|6|). Namely, G{x) is given by the statistical average Z' /Z, where the sum Z' 
is calculated in the same way as Z, but including the corresponding product of spin 
variables, which implies the following replacements: 



N 



(ri), 



(■ 



,(0 

x+l 



/ N 
\ i=i 



: case (a) (15) 
: case (b) . (16) 

The index i, entering in the sums Z' = J2i (r^^;^^ [case (a)] and Z' = J2i {j^2L+i) 
[case (b)], refers to the current configuration of the first row. These equations are 
obtained by an averaging in (|5|) and (P) over all the equivalent k values. Such a 
symmetrical form allows to reduce the amount of numerical calculations: due to the 
symmetry we need the summation over only ~ 2^ /N nonequivalent configurations 
of the first row instead of the total number of 2^ configurations. 

3.2 Improved algorithms 

The number of the required arithmetic operations can be further reduced if the 
recurrence relations (|8|) and ( |l3| ) are split into steps of adding single spin. To 
formulate this in a suitable way, let us first number all the 2-^ spin configurations 
{(t(1); ct(2); • • • ; a{N — 1); cr{N)} by an index i as follows: 



i = 1 

i =2 

i =3 

i =4 



{-1 
{-1 
{-1 
{-1 



-1; 
-1; 
-1; 
-1; 



{+1; +1; 



-1; -1; -1} 

-1; -1; +1} 

-1; +1; -1} 

; -1; +1; +1} 

; +1; +1; +1} 



(17) 
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Figure 2: Schematic pictures illustrating the algorithms of calculation for the lattices a and 
b introduced in Fig. ^ 

We remind that the sequence with k = 1,2, . . . , N corresponds to the numbers 

in the i-th row. The spin variables in (|l^) change just like the digits of subsequent 
integer numbers in the binary counting system. 

Consider now a lattice where n rows are completed, while the (n + l)~th row 
contains only £ spins where £ < N, as illustrated in Fig. |^ in both cases (a) and (b) 
taking as an example = 4. We consider the partial contribution (r„+i^). (i. e., 
i-th. component of vector r„,+i^£) in the partition function Z (or Z') provided by 
a fixed (i-th) configuration of the set of upper spins. These are the sequently 
numbered spins shown in Fig. |2| by empty circles. For simplicity, we have droped the 
index denoting the configuration of the first row. In case (b), the spin depicted by 
a double-circle has a fixed value a'. In general, this spin is the nearest bottom-left 
neighbour of the first spin in the upper row. According to this, one has to distinguish 
between odd and even n: a' refers either to the first (for odd n), or to the A'^-th 
(for even n) spin of the n-th row. It is supposed that the Boltzmann weights are 
included corresponding to the solid lines in Fig. |^ connecting the spins. In case (a) 
the weights responsible for the interaction between the upper numbered spins are 
not included. Obviously, for a given i > 1, Vn+i/ can be calculeted from rn+i/-i 
via summation over one spin variable, marked in Fig. |2| by a cross. In case (a) it 
is true also for £ = 1, whereas in case (b) this variable has fixed value a' at £ = 1. 
In the latter case the summation over a' is performed at the last step when the 
(n + l)-th row is already completed. These manipulations enable us to represent 
the recurrence relation (|8|) as 

Yn+l =Trn^WN Wn-1 • • • W^i r„ (18) 

with 

We= Weia) , (19) 

a=±l 

where the componets of the matrices Wi{a) are given by 

(Wiia))^^ = 5(j,ii(a,l,i))-exp(/3a{[cT(l)], + [a(2)], + [cT(Af)]J) 
(VF,(^))., = SU,h{a,£,i))-expi^a{[a{£)], + [ai£+l)]^}) : 1< ^ < A^ 
{WNia)),j = 6{j,ji{a,N,i))-exp{Pa[aiN)],) . (20) 
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Here S{j, k) is the Kronecker symbol and 

h{a,£,i) = i + ia-[ai£)]^)2''-'~' (21) 

are the indexes of the old configurations containing £—1 spins in the (n + l)-th row 
depending on the value a of the spin marked in Fig. ^a by a cross, as well as on the 
index i of the new configuration with £ spins in the (n + l)-th row, as consistent 
with the numbering ([T7|). 

The above equations (|l^) to (|20|) refers to ). In case (b) we have 

r„+i=ri,2r„^ ^ Wi^'^^wl^:^l--W^''^^wi''^\a')r^, (22) 

o-'=±l 

where ' are the matrices 

T^^(i'2) ^ ^ • (23) 

a=±l 

Here indexes 1 and 2 refer to odd and even row numbers n, respectively, and the 
components of the matrices ' (a) are 

(VFf •')(a))^^. = 5 {j,h,2{<^,£,i)) ■ exp {(3 [a{£)]^ {a + [a{e + 1)] J) , (24) 

where [ij{N + 1)]^ = a' and the index ji{a, £, i) is given by (^T|). For the other index 
we have 

J2(a,l,i) = 2i-2^-i([a(l)], + l) + i(a-l) 

i2(<T,Ai) = ji(fT,^,i) : ^>2. (25) 

— — ^(1 2) 

Note that the matrices We, and VF^ ' have only two nonzero elements in each 
row, so that the number of the arithmetic operations required for the construction 
of one row of spins via subsequent calculation of the vectors Tn+i,t increases like 
2N ■ 2^ instead of 2^^ operations necassary for a straightforward calculation of the 
vector rr„. Taking into account the above discussed symmetry of the first row, the 
computation time is proportional to 2^^L for both L x L [a] and \'^L x ^/2L (b) 
lattices in Fig. || with periodic boundary conditions. 

3.3 Application to different boundary conditions 

The developed algorithms can be easily extended to the lattices with antiperiodic 
boundary conditions. The latter implies that a{N + 1) = —a{N) holds for each 
row, and similar condition is true for each column. We can consider also the mixed 
boundary conditions: periodic along the horizontal axis and antiperiodic along the 
vertical one, or vice versa. To replace the periodic boundary conditions with the 
antiperiodic ones we need only to change the sign of the corresponding products of 
the spin variables on the boundaries. Consider, e. g., the case (a) in Fig. |. The 
change of the boundary conditions along the vertical axis means that the first term 
in the argument of the exponent in each of the Eqs. (|20| ) changes the sign for the last 
row, i. e., when n = L. The same along the horizontal axis implies that the term 
[(T(A^)]j in the equation for {Wi{a))^j changes the sign. In this case, however, the 
symmetry with respect to the configurations of the first row is partly broken and, 
therefore, we need summation over a larger number of nonequivalent configurations. 
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4 Transfer matrix study of critical Greens function and 
corrections to scaling in 2D Ising model 

4.1 General scaling arguments 

It is well known that in the thermodynamic limit the real-space Greens function of 
the Ising model behaves like G{r) oc r'^~'^~'^ at large distances r — > oo at the critical 
point f3 = [3c, where rj is the critical exponent having the value r/ = 1/4 in two 
dimensions (d = 2). Based on our transfer matrix algorithms developed in Sec ^, 
here we test the finite-size scaling and, particularly, the corrections to scaling at 
criticality. 

In the critical correlation function in the Fourier representation, i. e. G(k) at 
T = Tc, has been considered for the 99^ model. In this case the minimal value of the 
wave vector magnitude k is related to the linear system size L via kmin = 27r/L. In 
analogy to the consideration in Sec. 5.2 of p], on6 6xpGcts thctt k / kfji-iYi is an essential 
finite-size scaling argument, corresponding to r/L in the real space. In the Ising 
model at r ~ L one has to take into account also the anisotropy effects, so that the 
expected finite-size scaling relation for the real-space Greens function at the critical 
point j3 = (3c reads 

G{r) c^r'^-'^-'^ f{r/L) : r^oo,L^oo, (26) 

where the scaling function f{z) depends also on the crystallographic orientation of 
the line connecting the correlating spins, as well as on the orientation of the periodic 
boundaries. A natural extension of (pB]), including the corrections to scaling, is 

G{T) = Y,r-^^h{r/L), (27) 

where the term with Aq = d — 2 + is the leading one, whereas those with the 
subsequently increasing exponents Ai, A2, etc., represent the corrections to scaling. 
By a substitution fi{z) = z^' f'{z), the asymptotic expansion (^) transforms to 

G{r) = ar/L) L"^" + ^ L'^^ /,(r/L) j , (28) 

where fi{z) = f'f^{z) / fo{z) and uJi = \e — Aq are the correction-to-scaling exponents. 

We have tested the scaling relation ( [2^ ) in 2D Ising model by using the exact 
transfer matrix algorithms in Sec. ^. Consider, e. g., the correlation function in 
(10) crystallographic direction in the specific case (a) discussed in Sec. ^. According 
to (|2^), all points of f{r/L) = r^/^G(r) corresponding to large enough values of r 
and L should well fit a common smooth line, as it is actually observed in Fig ^ at 
2 < r < L/2 and L = 8, 12, 15, and 18. This example shows that the corrections 
to (^) are rather small, since the deviations from the spline curve constructed 
at L = 18 are very small at relatively small L values and even at r = 2. In 
general, our calculations provide a strong numerical evidence of correctness of the 
asymptotic expansion ([2^): at a given ratio r/L, the correlation function is described 
by an expansion in L powers with a really striking accuracy. Namely, our numerical 
analysis, which is sensitive to a variation of G{r) in the fourteens digit, does not 
reveal an inconsistency with this asymptotic expansion. 
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Figure 3: Test of the finite-size scaling relation G{r) ~ r^^^'^jirjV) in 2D Ising model on 
a square L x L lattice with periodic boundary conditions. The values of the scaling function 
f(r/L) estimated at L = 8 (triangles), L = 12 (circles), L = 15 (squares), and L = 18 
(pluses) well lie on a common smooth line showing that corrections to this scaling relation 
are small. 



4.2 Correction— to— scaling analysis for the L x L lattice 

Based on the scaling analysis in Sec. |4.1| , here we discuss the corrections to scaling 
for the lattice in Fig. |l^. We have calculated the correlation function G{r) at a fixed 
ratio r/L = 0.5 in (10) direction, as well as at r/L = 0.5\/2 in (11) direction at 
L = 2,4, 6, . . . with an aim to identify the correction exponents in ([2^). Note that 
in the latter case the replacement ([l^ ) is valid for G{\/2x) (where x = 1,2,3, .. .) 
with the only difference that A(x) = x. 

Let us define the effective correction-to-scaling exponent uJeff{L) in 2D Ising 
model via the solution of the equations 

L^/^G(r = const ■L) = a + h L'^^'^fi (29) 

at L = L, L + AL, L + 2AL with respect to three unknown quantities uJeff, cl, and 
b. According to (l28|), where Xq = r] = 1/4, such a definition gives us the leading 
correction-to-scaling exponent a; at L ^ oo, i. e., lim^^oo '^e//(i) = ^■ 

The calculated values of G{r = c - L) in the (10) and (11) crystallographic direc- 
tions [in case (a)] with c = 0.5 and c = 0.5\/2, respectively, and the corresponding 
effective exponents u;^ff(L), determined at AL = 2, are given in Tab. |l|. In both 
cases the effective exponent uj^ff{L) seems to converge to a value about 2. Besides, 
in the second case the behavior is more smooth, so that we can try someway to 
extrapolate the obtained sequence of Wg// values (column 5 in Tab. H) to L = oo. 
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Table 1: The correlation function G{r = c • L) in (10) (c = 0.5) and (11) (c = 0.5-\/2) crys- 
tallographic directions vs the linear size L of the lattice (a) in Fig. |l], and the corresponding 
effective exponents ujeffiL) and S3(L). 



L 


direction (10) 


direction (11) 


G(0.5L) 




G(0.5\/2L) 






2 


0.84852813742386 


2.7366493 


0.8 


1.8672201 




4 


0.74052044609665 


2.9569864 


0.71375464684015 


2.2148707 




6 


0.67202206468538 


1.8998036 


0.65238484475089 


2.1252078 




8 


0.62605120856389 


1.5758895 


0.60935351016910 


2.0611362 


1.909677 


10 


0.59238112628953 


1.6617494 


0.57724041054810 


2.0351831 


1.996735 


12 


0.56615525751968 


1.7774398 


0.55200680271678 


2.0232909 


2.002356 


14 


0.54485584658226 


1.8542943 


0.53141907668442 


2.0167606 


2.001630 


16 


0.52703456475995 




0.51414720882560 






18 


0.51178753041103 




0.49934511003360 







For this purpose we have considered the ratio of two subsequent increments in Wg/ / , 

0Jeff{L)-uJeff{L-AL) 
A simple analysis shows that r{L) behaves like 

r(L) = 1 - AL • (w' + 1)L~^ + o (i'^j (3^) 

at L — > 00 if u)eff{L) = to + (^L~^ ^ holds with an exponent to' > \. The numerical 
data in Tab. || show that Eq. (|3l|) represents a good approximation for the largest 
values of L at = 2. It suggests us that the leading and the subleading correction 
exponents in (|2^ ) could he uj = uJi = 2 and uj2 = 4, respectively. Note that LOeff{L) 
can be defined with a shift in the argument. Our specific choice ensures the best 
approximation by (|3^) at the actual finite L values. 

Let us now assume that the values of cje/ / {L) are known up to L = Lmax ■ Then 
we can calculate from (^) the r{L) values up to L = L^ax — and make a suitable 
ansatz like 

r{L) = 1 - SAL ■ L-'^ + b L-^ at L > Lmax (32) 



for a formal extrapolation of uJeff{L) to L = 00. This is consistent with ( [3l| ) where 
cj' = 2. The coefficient h is found by matching the result to the precisely calculated 
value at L = Lmax — AL. The subsequent values of LJeff{L), calculated from (30) 
and (32) at L > Lmax, converge to some value ui{Lmax) at L — > cxd. If the leading 
correction-to-scaling exponent uj is 2, indeed, then the extrapolation result u){Lmax) 
will tend to 2 at Lmax — > co irrespective to the precise value of uj' . 

As we see from Tab. the values of u){L) come remarkably closer to 2 as com- 
pared to iO(,ff{L), which seems to indicate that uj = 2. However, as we have discussed 
in Sec. there should be a nontrivial correction in (^) with uj = rj = 1/4. The 
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Table 2: The correlation function G{r — L) in (10) crystallographic direction and the 
effective exponents uJeff{L) and oj{L) vs the linear size L of the lattice (b) in Fig. 0. 



L 


G{L) 


We//(L) 




2 
3 


0.8 

0.7203484812087670 








4 


0.6690636562097066 








5 


0.6321925914229602 








6 


0.6037455936471098 








7 

8 


0.5807668304926868 
0.5616046762441826 


2.066235298 






9 


0.5452468033693456 


2.043461090 






10 


0.5310294874153481 


2.030235674 


1 


996772124 


11 


0.5184950262041604 


2.022130104 


1 


999333324 


12 


0.5073151480587211 


2.016864947 


1 


999941357 


13 


0.4972468711401118 


2.013265826 


2 


000036957 


14 


0.4881056192765374 


2.010701166 


2 


000040498 


15 


0.4797481011874659 


2.008811505 


2 


000044005 


16 


0.4720609977942179 


2.007380630 


2 


000053415 


17 


0.4649532511721054 


2.006272191 


2 


000063984 


18 


0.4583506666254706 


2.005396785 


2 


000073711 


19 


0.4521920457268738 








20 


0.4464263594840965 









fact that this nontrivial correction does not manifest itself in the above numerical 
analysis can be understood assuming that this correction term has a very small 
amplitude. 

4.3 Correction— to— scaling analysis for the ^/2L x \/2L lattice 

To test the existence of nontrivial corrections to scaling, as proposed at the end of 
Sec. |4.2| , here we make the analysis of the correlation function G{r) in (10) direction 
on the -v/2L x \/2L lattice shown in Fig. [l|b. The advantage of case (b) in Fig. |l] as 
compared to case (a) is that \/2 times larger lattice corresponds to the same number 
of the spins in one row. Besides, in this case we can use not only even, but all lattice 
sizes to evaluate the exponent to from calculations of G{r = L), which means that 
it is reasonable to use the step AL = 1 to evaluate tOeff and io{L) from Eqs. (p9D, 
(^) and (^. The results, are given in Tab. |. 

It is evident from Tab. |^ that the extrapolated values of the effective correction 
exponent, i. e. cj(-L), come surprisingly close to 2 at certain L values. Besides, 
the ratio of increments r [cf. Eq. (pOD] in this case is well approximated by (|3^), as 
consistent with existence of a correction term in ( p8|) with exponent 4. On the other 
hand, we can see from Tab. Q that Ati3(L) = u){L) — 2 tends to increase in magnitude 
at L > 13. We have illustrated this systematic and smooth deviation in Fig. ^ The 
only reasonable explanation of this behavior is that the expansion ( |28| ) necessarily 



12 
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4e-05 
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-8e-05 ' — — — — ' — — — — — — — . 

0.0 0.2 0.4 (10/Lj 

Figure 4: The deviation of the extrapolated effective exponent Auj{L) = — 2 as a 
function of L~*. The extrapolation has been made by using the calculated G(r) values in 
Tab. ^ up to the size L + 2. A linear convergence to zero would be expected in absence of 
any correction term with exponent lu < 2. 



contains the exponent 2 and, likely, also the exponent 4, and at the same time it 
contains also a correction of a very small amplitude with u < 2. The latter explains 
the increase of Aa;(L). Namely, the correction to scaling for L^/^G{L) behaves like 
const ■ L^^ [l + o (L~^) +eLp^~^^^ with e ^ 1, which implies a slow crossover of 
the effective exponent LO(.ff{L) from the values about 2 to the asymptotic value w. 
Besides, in the region where eL^~^ <^ 1 holds, the effective exponent behaves like 

We/ / (L) ~ 2 + biL^-"^ + , (33) 

where 61 <C 1 and 62 are constants. By using the extrapolation of We// with u)' = 2 
in (|3l|) and (|3^), we have compensated the effect of the correction term 62^^^. Be- 



sides, by matching the amplitude h in (32) we have compensated also the next trivial 
correction term ~ in the expansion of uJeff{L). It means that the extrapolated 
exponent u}{L) does not contain these expansion terms, i. e., we have 

w(L) = 2 + biL^-'^ + 5u{L) , (34) 

where 5u){L) represents a remainder term. It includes the trivial corrections like 
L~^, L~^, etc., and also subleading nontrivial corrections, as well as corrections of 
order (eL^""^)^, [eL'^~'^)^, etc., neglected in (|3^). According to the latter, Eq. (|3^ 
is meaningless in the thermodynamic limit L — > 00, but it can be used to evaluate 
the correction-to-scaling exponent lo from the transient behavior at large, but not 
too large values of L where hiLp'^^ <^ 1 holds. In our example the latter condition 
is well satisfied, indeed. 
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0.00 



0.08 



(10/L) 



5.75 



0.16 



Figure 5: The exponent 2 — uj estimated from (^7|) at different system sizes. From top 
to bottom (if looking on the left hand side): a = 0,1,3.5,5.75,7.243,9.25,12. The resuhs 
at the optimal a value 7.243 are shown by empty circles. The dashed line indicates our 
theoretical asymptotic value 2 — lu ~ 1.75, whereas the dot-dashed line - that proposed 
in 



Based on (^4|), we have estimated the nontrivial correction-to-scaling exponent 
(jj by using the data of cj(L) in Tab 0. We have used two different ansatzs 



2-cJi(L) =ln [Au{L)/Ai:!{L- l)]/ln[L/(L- 1)] 



and 



2 - uJ2{L) = L [Aw(L) - Ai:i{L - 1)] /At^(L) , 
as well as the linear combination of them 

uj{L) = (1 — a) uJi{L) + a uJ2{L) 



(35) 
(36) 

(37) 



containing a free parameter a. We have uj{L) = uJi{L) at a = and io{L) = uj2{L) 
at a = 1. In general, the effective exponent td{L) converges to the same result w at 
arbitrary value of a, but at some values the convergence is better. The results for 
2 — (jj{L) vs L^~^ at different a values are represented in Fig. ^ by a set of curves. In 
this scale the convergence to the asymptotic value would be linear (within the actual 
region where L 3> 1 and bi <C 1 hold) for a = at the condition 5(jj{L) oc L~^. 
We have choosen the scale of L~^'^^, as it is consistent with our theoretical prediction 
in Sec. ^ that lo = 1/4. Nothing is changed essentially if we use slightly different 
scale as, e. g., L"^"^/^ consistent with the correction-to-scaling exponent u = 4/3 
proposed in |17]. As we see from Fig. ^, all curves tend to merge at our asymptotic 
value 2 — IV = 1.75 shown by a dashed line. The optimal value of a is defined by 
the condition that the last two estimates w(17) and w(18) agree with each other. It 
occurs at a = 7.243, and the last two points lie just on our theoretical line. 
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It is interesting to compare our results with those of the high temperature (HT) 



series analysis in The authors of have found "almost by inspection" a 



correction with exponent u; = 9/4 in the asymptotic expansion of the susceptibility 
X- If such a correction exists in the susceptibility, it must be present also in the cor- 
relation function due to the relation x = X]x^(^i ~^)- Surprisingly, our extremely 
accurate calculations by exact algorithms have not revealed such a correction. Our 
analysis shows that nontrivial corrections exist, indeed, in the correlation function 
of 2D Ising model, but they are so extremely small and strongly masked by triv- 
ial corrections like and that they could not be detected by approximate 
methods like HT series expansion. 

4.4 Comparison to the known exact results and estimation of nu- 
merical errors 

We have carefully checked our algorithms comparing the results with those obtained 
via a straightforward counting of all spin configurations for small lattices, as well as 
comparing the obtained values of the partition function to those calculated from the 
known exact analytical expressions. Namely, an exact expression for the partition 
function of a finite-size 2D lattice on a torus with arbitrary coupling constants 
between each pair of neighbouring spins has been reported in jl^ obtained by the 
loop counting method and represented by determinants of certain transfer matrices. 
In the standard 2D Ising model with only one common coupling constant (3 these 
matrices can be diagonalized easily, using the standard techniques |^0|. Besides, the 
loop counting method can be trivially extended to the cases with antiperiodic or 



mixed boundary conditions discussed in Sec. 3.3. It is necessary only to mention 
that each loop gets an additional factor —1 when it winds round the torus with 
antiperiodic boundary conditions. We consider the partition functions Zpp = Z, 
Zaa^ Zap, Zpa- In this notation the first index refers to the horizontal or x axis, 
and the second one - to the vertical or y axis of a lattice illustrated in Fig. |l]a; p 
means periodic and a ~ antiperiodic boundary conditions. As explained above, the 
standard methods leads to the following exact expressions: 

Zpp = (Qi + Q2 + Q3-Qo)/2 

Zap = (Qo + Ql+Q3-Q2)/2 

Zpa = {Qo + Ql+Q2-Q3)/2 (38) 
Zaa = {Qo + Q2 + Q3-Ql)/2 

where Qq is the partition function represented by the sum of the closed loops on the 
lattice, as consistent with the loop counting method in ||2^, whereas Qi, Q2, and Q3 
are modified sums with additional factors exp(Ax-i7r/A^+ Ay-ivr/L), exp(Ax-i7r/A^), 
and exp(Ay •■ivr/L), respectively, related to each change of coordinate x by Ax = ±1, 
or coordinate y by Ay = ±1 when making a loop. The standard manipulations ||20[| 
yield 

Q. = 2^^ Yl [cosh2(2/3) - sinh(2/3) (39) 

1/2 



X cos 



cos 
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where the wave vectors qx = {2tt /N) ■ n and qy = (27r/L) • i run over all the values 
corresponding to n = 0, 1, 2, . . . , — 1 and £ = 0, 1, 2, . . . , L — 1. In the case of 
the periodic boundary conditions, each loop of Qq has the sign \19\^ 



where m is the number of intersections, a is the number of windings around the 
torus in x direction, and 6 - in y direction. The correct result for Zpp is obtained if 
each of the loops has the sign (—1)"^. In all other cases, similar relations are found 
easily, taking into account the above defined additional factors. Eqs. (|3^) are then 
obtained by finding such a linear combination of quantities Qi which ensures the 
correct weight for each kind of loops. 

All our tests provided a perfect agreement between the obtained values of the 
Greens functions G{r) (a comparison between straightforward calculations and our 
algorithms), as well as between partition functions for different boundary conditions 
(a comparison between our algorithms and Eq. (|38|)). The relative discrepancies were 
extremely small (e. g., 10~^^), obviously, due to the purely numerical inaccuracy. 

We have used the double-precision FORTRAN programs. The main source of 
the inaccuracy in our calculations is the accumulation of numerical errors during 
the summation of long sequences of numbers, i. e., during the sumation over all the 
nonequivalent configurations of the first row of spins. To eliminate the error for the 
largest lattice L = 20, we have split the summation in several parts in such a way 
that a relatively small part, including only the first 10 000 configurations from the 
total number of 52 487 nonequivalent ones, gives the main contribution to Z and 
Z' . The same trick with splitting in two approximately equal parts has been used 
at L = 19, as well. By comparing the summation results with different splitings, 
we have concluded that a systematical error in G{r) at L = 20 could reach the 
value about 3 • 10^^^. The calculations at L = 18 and 19 have been performed with 
approximately the same accuracy. The systematical errors in subsequent G{r) values 
tend to compensate in the final result for lo{L). The resulting numerical errors in 
Fig. ^ are about 0.02 or 0.03, i. e., approximately within the symbol size. 

5 Analysis of the partition function zeros in 3D Ising 
model 

In this section we discuss the recent MC results |^] for the complex zeros of the 
partition function of the three-dimensional Ising model. Namely, if the coupling (5 
is a complex number, then the partition sum has zeros at certain complex values of 
/3 or u = . The nearest to the real positive axis values j3i and are of special 
interest. Neglecting the second-order corrections, ti? behaves like 



+ AL-i/'^ + BL-(i/'^)-^ (40) 



at large L, where Uc = e is the critical value of u, A and B are complex constants, 
and Lo is the correction-to-scaling exponent. According to the known solution given 



in [ pO| , 21 1, the partition function zeros correspond to complex values of sinh(2/3) 
located on a unit circle in the case of 2D Ising model, so that A is purely imaginar. 
This solution, however, is only asymptotically exact at L ^ oo. Nevertheless, based 
on an analysis of the exact expression Eq. (|3^) we conclude that the statement 
Re A = is correct. This fact is obvious in the case of Brascamp-Kunz boundary 
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(1 /\^+(0 (RG value) 




Figure 6: Effective critical exponents for the real (empty circles) and the imaginary (solid 
circles) part of complex partition-function-zeros of 3D Ising model depending on L^^/^, 
where L is the linear size of the system. Solid lines show the linear least-squares fits. The 
asymptotic values from our theory are indicated by horizontal dashed lines, whereas those 
of the RG theory - by dot-dot-dashed lines. A selfconsistent extrapolation within the RG 
theory corresponds to the tiny dashed line. 



conditions [p2| , p^ j. The latter means that the critical behavior of real and imaginary 
parts of u\ — Uc essentially differ from each other, i. e.. Re (ui — Uc) oc L~^^^'^^^^ and 
Im [ui) oc L^^/'^ (where, in this case of d = 2, u = uj = 1) at L ^ oo. The MC data 
of in fact, provide a good evidence that the same is true in three dimensions. 

Based on MC data for the partition function zeros in 3D Ising model, the authors 
of Ref. have searched the way how to confirm the already known estimates for 
v. Their treatment, however, is rather doubtful. First, let us mention that, in 
contradiction to the definition in the paper, Ui values listed in Tab. I of |^] are 
not equal to e~^^ (they look like e~^^i). Second, the fit to a theoretical ansatz for 
I u\{L) — Uc\, Eq. (6) in |^], is unsatisfactory. This ansatz contains a mysterious 
parameter 03. If we compare Eqs. (5) and (6) in Q, then we see immediately that 
03 = {\/v) + uj. At the same time, the obtained estimate for 03, i. e. 03 = 4.861(84), 
is completely inconsistent with the values of (l/i^) + w, about 2.34, which follow 
from authors own considerations. Our prediction, consistent with the correction- 
to-scaling analysis in Sec. § (and with £ = 4 in (^) to coincide with the known exact 
result at d = 2), is = 2/3 and u) = 1/2, i. e., {l/p) + to = 2. 

To obtain a more complete picture, we have considered separately the real part 
and the imaginary part of u\ — Uc- We have calculated from f3^ data listed in 
Tab. I of 1^ and have estimated the effective critical exponents y'^jj(L) and y'ljj{L), 
separately for Re {u\ — u^) and Im (n^) , by fitting these quantities to an ansatz 
const ■ L~^'^ff and const • L'^^ff, respectively, at sizes L and L/2. The value of 
Uc consitent with the estimation of the critical coupling in [^], f3c — 0.2216545, 
has been used. The results are shown in Fig. As we see, y^jj (empty circles) 
claims to increase above y'^jjr (solid circles) when L increases. This is a good nu- 
merical evidence that, like in the two-dimensional case, the asymptotic values are 
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y' = limL^ooye//(^) = (1/^) ^'"d y" = liniL^oo yg// (-^) = 1/^- According to 
our theory, the actual plots in the L~^/^ scale are linear at L ^ oo, as consistent 
with the expansion in terms of L'^^ . The linear least-squares fits are shown by solid 
lines. The zero intercepts 1.552 and 1.899 are in approximate agreement with our 
theoretical values 1.5 and 2 indicated by horizontal dashed lines. The relatively small 
discrepancy, presumably, is due to the extrapolation errors and inaccuracy in the 
simulated data. The result for y' is affected by the error in /3c value. However, this 
effect is negligibly small. Assuming a less accurate value /3c = 0.221659, consistent 



with the estimations in |32, 30|, we obtain y' = 1.914. 

The behavior of y'^jrj is rather inconsistent with the RG predictions. On the one 
hand, y'^jjr claims to increase above y'^jjr and also well above the RG value of l/z^ 
(the lower dod-dot-dashed line at 1.586), and, on the other hand, the extrapolation 
yields y' value (1.899) which is remarkably smaller than (1/z^) + ~ 2.385 (the 
upper dot-dot-dashed line) predicted by the RG theory. For selfconsistency, we 
should use the linear extrapolation in the scale of L~'^ with tu = 0.799 (the RG 
value). However, this extrapolation (tiny dashed line in Fig. |6|), yielding y' ~ 1.747, 
does not solve the problem in favour of the RG theory. 

The data points of y^^j look (and are expected to be) less accurate than those of 
y'effj since Re (m? — Uc) has a very small value. The y'^jij data do not look scattered, 
therefore they allow a refined analysis with account for nonlinear corrections. To 
obtain stable results, we have included the data for smaller lattice sizes L = 3 and 
L = 4 given in p3]. In principle, we can use rather arbitrary analytical function 



0(/3) to evaluate the effective critical exponent 

y'^jfiL) = In [/mc/. (/30(L/2)) /Im (/3?(L))] / In 2 

and estimate its asymptotic value y" . For an optimal choice, however, y"jj-{L) vs 
L~'^ plot should be as far as possible linear to minimize the extrapolation error. In 



this aspect, our choice (j) = exp(— /3) is preferable to (f) = exp(— 4/?) used in |24|| . 
We have tested also another possibility, i. e. (p = sinh(2/?), which appears as a 
natural parameter in the case of 2D Ising model. The shape of the y"jj-{L) plot 
can be satisfactory well approximated by a third-order, but not by a second-order, 
polinomial in as it can be well seen when analyzing the local slope of this 

curve. The corresponding four parameter least-squares fits are shown in Fig. 0. They 
yield y" ~ 1.473 in the case of cj) = exp(— /3) (long-dashed line) and y" ~ 1.518 at (j) = 
sinh(2/3) (solid line). It is evident from Fig. |^that in the latter case we have slightly 
better linearity of the fit, therefore l/i^ ~ 1.518 is our best estimate of the critical 
exponent 1 /z/ from the actual MC data. Thus, while the row estimation provided the 
value y" = 1/z/ ~ 1.552 which is closer to the RG prediction 1/zv ~ 1.586 (horizontal 
dot-dot-dashed line), the refined analysis reveals remarkably better agreement with 
our (exact) value l/v = 1.5 (horizontal dashed line). 
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Figure 7: Effective critical exponent y'^j.j{L) for the imaginary part of the complex 
partition-function-zeros as a function of L^^/^, where L is the linear size of the system. The 
empty circles correspond to cf) = exp(— whereas the solid circles to (f) = sinh(2/3). The 
corresponding least-squares fits y'^ff{L) = 1.4734 + 1.3321L-1/2 - 4.7587L-1 + 6.8894i-3/2 
and y'effiL) = 1.5184 + 0.7271i-i/2 _ 2.0309i-i + 3.3095X^3/2 ^^.^ gj^^^^^ long-dashed 
line and solid line, respectively. Our asymptotic value y" = l/i/ = 1.5 is indicated by 
horizontal dashed line, whereas that of the RG theory (1.586) - by dot-dot-dashed line. 



6 Xcp'^ model and its crossover to Ising model 

Here we discuss a (/?^ model on a three-dimensional cubic lattice. The Hamiltonian 
of this model, further called Xip'^ model, is given by 

H/T = Y,l-2Kj2ip^^^+f, + ^l + x[^l-iy\ , (41) 

where the summation runs over all lattice sites, T is the temperature, 
(^x £ ] — cxd;+oo[ is the scalar order parameter at the site with coordinate x, ft 
is a unit vector in the //-th direction, k and A are coupling constants. Obviously, 
the standard 3D Ising model is recovered in the limit A — > oo where fluctuations 
are suppressed so that, for a relevant configuration, cp'^ I oi ip^ ^ ±1 holds. The 
MC data for the Binder cumulant in this Xip^ model have been interpreted in ac- 
cordance with the e-expansion and a perfect agreement with the conventional RG 
values of critical exponents has been reported in [25|. According to the definition 



in 1 25 1, the Binder cumulant U is given by 



{m? 



f^Si. (12) 



where m = L ^ (^x is the magnetization and L is the linear size of the system. 



Based on the e~expansion, it has been suggested in |25| that, in the thermodynamic 
limit L — > oo, the value of the Binder cumulant at the critical point n = Kc(A) 
and, equally, at a fixed ratio Za/Zp = 0.5425 (the precise value is not important) of 
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partition functions with periodic and antiperiodic boundary conditions is a universal 
constant U* independent on A. We suppose that the latter statement is true, but 
not due to the e-expansion. It is a consequence of some general argument of the 
RG theory: on the one hand, U is invariant under the RG transformation and, 
on the other hand, an unique fixed point (not necessarily the Wilson-Fisher fixed 
point) exists in the case of an infinite system, so that U = U* holds at L — > c« and 
(A) where U* is the fixed-point value of U. The above conclusion remains true 
if we allow that the fixed point is defined not uniquely in the sense that it contains 
some irrelevant degree(s) of freedom (like c* and A in the perturbative RG theory 
discussed in Sec. 2 of ||5|) not changing U. The numerical results in confirm the 
idea that liniL^ao U (L) = U* holds at criticality, where U* is a universal constant 
independent on the specific microscopic structure of the Hamiltonian. 

7 Estimation of the correction exponent lu 

Based on the idea that U* is constant for a given universality class, here we estimate 
the correction-to-scaling exponent uj. According to Sec. ||, corrections to finite-size 
scaling for the magnetization of the actual 3D Ising and A(/?^ models are represented 
by an expansion in terms of L~'^ where = 1/2. One expects that the magnetization 
(Binder) cumulant (^) has the same singular structure. Since lim^^oo U{L, A) = U* 
holds at a fixed ratio Za/Zp, a suitable ansatz for estimation of uj is p5| 

U{L, Xi)-U{L,X2)^ const ■L-'' at Z^/Zp = 0.5425 , (43) 

which is valid for any two different nonzero values Ai and A2 of the coupling constant 
A. The data for AU{L) = U{L, 0.8) - U{L, 1.5) can be read from Fig. 1 in [|| (after 
a proper magnification) without an essential loss of the numerical accuracy, i. e., 
within the shown error bars. Doing so, we have evaluated the effective exponent 

cOeffiL) = In [AC/(L/2)/AC/(L)] / In 2 , (44) 

i. e., uJeff{l2) ~ 0.899, uj^ffilG) ^ 0.855, and u;e//(24) ~ 0.775. These values 
are shown in Fig. ^ by crosses. Such an estimation, however, can be remarkably 
influenced by the random scattering of the simulated data points, particularly, at 
larger sizes where AU{L) becomes small. This effect can be diminished if the values 
of AU{L) are read from a suitable smoothened curve. We have found that AU{L) 
within L G [7; 24] can be well approximated by a second-order polinomial in 
The corresponding refined values a;e//(16) ~ 0.8573 and a;ejj(24) ~ 0.7956 read 
from this curve are depicted in Fig. ^ by empty circles. These values are similar to 
those obtained by a direct calculation from the original data points (crosses). 

In such a way, we see from Fig. ^ that the effective exponent uJeffiL) decreases 
remarkably with increasing of L. According to GFD theory, iVeffiL) is a linear 
function of L^^/^ at L — > 00, as consistent with the expansion in terms of L^'^ 
where uj = 0.5. More data points, including larger sizes L, are necessary for a 
reliable estimation of the asymptotic exponent uj = lim^^oo '^e//(-^)- Nevertheless, 
already a row linear extrapolation in the scale of with the existing data points 

yields the result uj ~ 0.547 which is reasonably close to the exact value 0.5 (horizontal 
dashed line in Fig. |8|) found within the GFD theory. The corresponding least-squares 
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Figure 8: Effective correction-to-scaling exponent uJeff{L) in the 0(n)-synimetric X(p'^ 
model with n — I (empty circles and crosses) and 0(2)-symmetric dd — XY model (solid 
circles) depending on the system size L. The linear least-squares fits give row estimates of 
the asymptotic uj values 0.547 and 0.573, at n = 1 and 2, respectively. The corresponding 
theoretical values of the GFD theory 1/2 and 5/9 (used in the L~'^ scale of the horizontal 
axis) are indicated by arrows. The dot-dot-dashed line shows the value 0.845(10) proposed 
in for the 3D Ising universality class (n = 1). 



fit with circles (at L = 24, 16) and cross (at L = 12) is shown in Fig. ^ by a straight 
solid line. It is evident from Fig. ^ that the final result u = 0.845(10) (horizontal 
dot-dot-dashed line) reported in |25| represents some average effective exponent for 



the interval L S [6; 24]. It has been claimed in |2^ that the estimates for lo (cf. 
Tab. 2 in |^5|) are rather stable with respect to a variation of Lmim where Lmin is 
the minimal lattice size used in the fit. Unfortunately, the analysis has been made 
in an obscure fashion, i. e., giving no original data, so that we cannot check the 



correctness of this claim. Besides, the estimates in Tab. 2 of [25| has been made by 
using an ansatz 

U{L, \) = U* + ci(A)L-'^ at Za/Zp = 0.5425 , (45) 



which is worse than (|4^). Namely, ( |4^ ) and (45) are approximations of the same 



order, but (45) contains an additional parameter U* which is not known precisely. 
The results of an analysis with the ansatz (|4^), reflected in Tab. 5 of [^5[, are 
not convincing, since only very small values of Lmin (up to Lmin = 6) have been 
considered. 

In any case, we prefer to rely on that information we can check, and it shows that 
the claim in that ui = 0.845(10) holds with ±0.01 accuracy cannot be correct, 
since uJeff{L) is varied in the first decimal place. 

We have made a similar estimation of to for the dynamically diluted 0(2)- 
symmetric (n = 2) XY (dd—XY) model simulated in (n = 2). In the case of the 



dd—XY model, parameter D (cf. Eq.(6) in p6|]) plays the role of A in (|^). The data 
for the Binder cumulant in Fig. 1 of [p6| look rather accurate, i. e., not scattered. 
This enables us to estimate uJeff just from the data at D = 1.03 and D = oo (XY 
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model). The resulting values of ujeff are depicted in Fig. ^ by solid circles. The scale 
of L~'^ is used, where lo = 5/9 is our theoretical value of the correction-to-scaling 
exponent at n = 2 consistent with the general hypothesis proposed in Sec. ||. As 
we see, the solid circles can be well located on a smooth line which, however, is 
remarkably curved at smaller sizes. Due to the latter reason, we have used only 
the last three points (the largest sizes) for the linear fit (solid line) resulting in an 
estimate uj 0.573 which comes close to our theoretical value a; = 5/9 = 0.555 . . . 

In summary, the extrapolated u! values (Fig. ^ in both cases n = 1 and 2 
are reasonably close to our theoretical values 1/2 and 5/9 indicated by arrows. 
Only a small systematic deviation is observed. This, likely, is due to the error of 
linear extrapolation: the uJeff{L) plots have a tendency to curve down slightly. The 
conventional (RG) estimate uj ~ 0.8 more or less corresponds to effective exponents 
for currently simulated finite system sizes, but not to the asymptotic exponents. 



8 Fitting the susceptibility data at criticality 

In this section we discuss some fits of MC data at criticality. According to the 
finite-size scaling theory, the susceptibility x near the critical point is represented 
by an expansion 

X = L^-" (^go{L/0 + E L-'^'miL/O^ , (46) 

where gi{L/^) are the scaling functions, ^ is the correlation length of an infinite 
system, rj is the critical exponent related to the k~'^~^^ divergence of the correlation 
function in the wave vector space at criticality, and ivi are correction-to-scaling 
exponents, uji = uj being the leading correction exponent. The correlation length 
diverges like ^ oc t~'^ at t — > 0, where t = 1 — /5//3c is the reduced temperature. Thus, 
for large L, in close vicinity of the critical point where tL^^'^ ^ 1 holds Eq. ( |4^ ) can 
be written as 

X = a L^-^ + E hL-'^' + Kt. L)^ , (47) 

where a = go(0) and bi = gi{0)/go{0) are the amplitudes, and 6(t,L) is a correction 
term which takes into account the deviation from criticality. In the first approxima- 
tion it reads 

6{t, L) ^ c ■ tL^/" , (48) 

where c is a constant. 

We start our analysis with the standard 3D Ising model with the Hamiltonian 

H/T = -PY.aiaj . (49) 

(ij) 



The critical point of this model has been found in [27| to be (ic — 0.2216545. We have 



made our own tests with the data of [27|, and have obtained the same value within 
the uncertainty of ±10"^. From the maximal values of the derivative dhi{m?) / df5 = 
dhix/d(3 evaluated in |^0| we conclude that the shift of [3 by 10~^ produces the 



variation of In x at L = 96 near (3 = j3ci which does not exceed 4.7-10 in magnitude 
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Figure 9: The fits of In (x/^^) data at criticality (ansatz (|47|)) shifted by a constant c. SoHd 
circles represent the MC data for 3D Ising model |2|l at /3 = 0.2216545 (c = 0). The fits 
with our (GFD) exponents (Ina = 1.065289, 6i = -2.72056,^2 = 8.18686,63 = -10.49614) 
and with those of |||2| (Ina = 0.430933, 61 = 0.05850,^2 = -7.74767,63 = 12.42890) 
are shown by solid and tiny-dashed lines, respectively. The empty boxes are MC data for 
3-component 3D XY model pi, shifted by c = 0.85. 



The latter means that, with a good enough accuracy, we may assume that (3c is just 
0.2216545 when fitting the susceptibility data at criticality within L G [4; 128]. Here 
we mean the MC data given in Tab. 25 of |^^. We have made and compared several 
fits of these data to ansatz (^) with 6{t, L) = (more precisely, to the corresponding 
formula for In^) for two different sets of the critical exponents, i. e., our (GFD) and 
that proposed in [^]. The fits made with our exponents systematically improve 
relative to those made with the exponents of [^], as the system sizes grow and the 
approximation order increases. The necessity to include several correction terms is 
dictated by the fact that corrections to scaling are rather strong. According to the 
least-squares criterion, the fit with our exponents r] = 1/8 and uji = 1/2 becomes 
better than that provided by the more conventional exponents rj = 0.0358(4), uji = 
0.845(1), UJ2 = 2wi, and ^3 = 2 |^] starting with Lmin = 28 (i. e., L G [Lmin, 128]), 
if two correction terms (/ = 1, 2) are included. In the case of three correction terms 
it occurs already at Lmin = 11- The four-parameter (a, bi, 62, 63) fits to MC data 
(empty circles) within L S [14; 128] are shown in Fig. |^. The fit with our exponents 
(upper solid line) is relatively better at larger sizes. However, both fits (upper solid 
and dashed lines) look, in fact, quite similar, so that we cannot make unambiguous 
conclusions herefrom. 

For comparison, we have shown in Fig. ^ also the MC data for 3D XY model |29], 
where the order parameter is 3-component vector with only two interacting compo- 
nents. As we see, the actual MC data (empty boxes) at Pc evaluated approximately 
Pc — 0.6444 |29| are rather scattered and, therefore, unsuitable for a refined analy- 
sis. Nevertheless, this is a typical situation where authors of such data make a very 
"accurate" and "convincing" estimation 7/1/ = 1.9696(37) or rj = 0.0304(37) making 
a simple linear fit. However, the refined analysis given above has shown that even 
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Figure 10: The effective critical exponent rieff{L) (solid circles) obtained by fitting the 
susceptibility data of 3D Ising model at criticality (/3 — 0.2216545) |^ within the inter- 
val [L;2L]. The least-squares approximations obtained by fitting the rjeff{L) data within 
[Lmi„;64] to a third-order polinomial in are shown by dashed {Lmin = 9), solid 

{Lmin = 10), and tiny-dashed {Lmin — 12) lines. The asymptotic value rj = 1/8 of the GFD 
theory is indicated by a horizontal dashed line. The dot-dot-dashed line represents the rj 
value 0.0358 proposed in [psf . 



in the case of 3D Ising model, where the data are incompatibly more accurate, it 
is not so easy to distinguish between r] = 0.0358 and r] = 1/8. Moreover, a refined 
analysis prefer the second value which is much larger than those usually provided by 
linear fits at typical system syzes L < 48. This is particularly well seen in Fig. |l^, 
where the effective critical exponent r](,ff{L) of the 3D Ising model, estimated via 
the linear fit within [L; 2L], is depicted by solid circles. As we see, r]eff{L) tends to 
increase well above the conventional value 0.0358 (horizontal dot-dot-dashed line). 
The shape of the rieff{L) plot is satisfactory well reproduced by a third-order poli- 
nomial in the actual scale of L~^/^. Three such kind of least-squares approximations 
(at Lmin = 9, 10, 12) are shown in Fig. IC. These fits do not provide very accurate 
and stable asymptotic values of ij. Nevertheless, they are more or less in agreement 
with our theoretical prediction r] = 1/8 (horizontal dashed line). Besides, the values 
of 7]eff are affected by the error in /3c (about 10"'') only slightly, i. e., by an amount 
not exceeding 0.001. 



9 A test for 3D Ising model with "improved" action 

Here we discuss some estimations of the critical exponents from the susceptibility 
data of 3D Ising model, reported in p^], with the so called "improved" action (i. e., 
LL/T). One of the problems with the standard 3D Ising model is that corrections to 
scaling are strong. It has been proposed in |27| to solve this problem by considering 
a modified (spin-1) Ising model with the Hamiltonian 

H/T = -/3 E ^^^i + ^ E ' (50) 
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Figure 11: The effective amplitudes 106i (circles) and 62 (rhombs) in ( ^ ) estimated at 
fixed exponents rj = 0.0358, uji = 0.845, lu2 = 2uji, and 1^ = 0.6296 by fitting the MC 
data within L S [Lmin',56]. Filled symbols correspond to S{t,L) = 0, empty symbols - to 
(5(t, L) = 10~^L^/'^. The effective amplitudes bi and &2 estimated with the critical exponents 
of our GFD theory (77 = 1/8, = 1/2) at S{t, L) = are shown by "x" and "+" , respectively. 
Lines represent the least-squares approximations by a fourth-order polinomial in L. 



where the spin cTj takes the values 0, ±1, with two coupling constants /3 and D ad- 
justed in such a way that the leading correction to finite-size scaling vanishes for all 
relevant physical quantities (magnetization cumulant, energy per site, susceptibility, 
etc.) and their derivatives. Moreover, according to the claims in [p7| (see the con- 
clusions in |2^), the ratios of the leading and subleading corrections are universal, 
so that not only the leading but all (!) corrections should vanish simultaneously. 

We have checked the correctness of these claims as described below. We have 
fitted the corresponding to ( |47|) expression for In x to the susceptibility data of the 
"improved" 3D Ising model (|g) with (/3,L>) = (0.383245,0.624235) (this is an ap- 
proximation of the critical point) given in ||2^ (Tab. 26). By fixing the exponents, 
the least-squares fit within L G [Lmin, 56] (here L = 56 is the maximal size available 
in Tab. 26 of ||2^), including the leading and the subleading correction to scaling, 
provides the effective amplitudes a, 61, and 62 depending on Lmin- We have made 
a test with the critical exponents rj = 0.0358(4), uj = 0.845(10), and 1^ = 0.6296(3) 
proposed in [p^]. These values are close to those of the usual RG expansions ||28[|, 



but, as claimed in |25], they are more accurate. According to [25], the asymptotic 
expansion contains corrections like L"""^ and L^^", where n = 1,2,3,... Thus we 
have u>i = uj and lj2 = 2uj. The resulting amplitudes Wbi{Lmin) and h2{Lmin) are 
shown in Fig. |ll|by circles and rhombs, respectively. We have depicted by filled sym- 
bols the results of the fitting with 5(t, L) = 0, assuming that the critical coupling 
(ic = 0.383245 has been estimated in [^5| with a high enough (6 digit) accuracy. The 
data points quite well fit smooth (tiny dashed) lines within Lmin G [4; 20] , which 
means that the statistical errors are reasonably small. If the exponents used in the 
fit are correct and corrections to scaling are small indeed, then the convergence of 
the effective amplitudes to some small values is expected with increasing of Lmin- 



25 



However, as we see from Fig. 11, the effective amplitudes tend to increase in mag- 
nitude acceleratedly as Lmin exceeds 14. A small inaccuracy in (3c value can be 
compensated by the term 5{t, L) ~ c*L}/^ in (|^), where c* = ct (cf. Eq. (|48|)). The 
results of fitting with c* = 10~^ are shown in Fig. |ll| by empty symbols. As we see, 
the expected inaccuracy in /3c of order 10~^ does not change the qualitative picture. 
The increase of the effective amplitudes indicates that either the exponents are false, 
or the asymptotic amplitudes are not small (or both). This is our argument that the 
claims in about very accurate critical exponents, extracted from the 3D Ising 
model with "improved" action, are incorrect. 

For comparison, we have shown in Fig. |11| also the effective amplitudes ^i(-^rn.m) 
and b2{Lmin) (by "x" and "+", respectively) estimated with the critical exponents 
of our GFD theory (Sec. D (r/ = 1/8, uji = 1/2), assuming 6{t, L) = 0. The effective 
amplitudes converge to some values with increasing of Lmin- These, however, are not 
the true asymptotic values, since the maximal size of the system has been eliminated 
to L = 56. 



10 A test for the standard 3D Ising model 

A test with the effective amplitudes, as in Sec. appears to be more sensitive tool 
as compared to the fits discussed in Sec. @. Here we consider the standard 3D Ising 



model. We have fitted all data points in Tab. 25 of [^7[ within the interval of sizes 
[L; 8L] to the theoretical expression for In^ (consistent with ([47|)) to evaluate the 
effective amplitudes a and bi with / = 1,2,3 depending on L. Exceptionally in 
the case if all the involved exponents are correct (exact) each effective amplitude 
can converge to a certain nonzero asymptotic value at L — > oo. In other words, if 
one tries to compensate the inconsistency in the exponent by choosing appropriate 
amplitude, then the amplitude tends either to zero or infinity at L — > cxd. 

We have shown in Fig. |l^ the effective amplitudes Ina(-L) and bi{L) in the case 
of our critical exponents r/ = 1/8 and oji = 1/2. As we expected, the effective 
amplitudes converge to some nonzero values with increasing of L. This is a good 
numerical evidence that our critical exponets are true. The case with the exponents 
of H] ?7 = 0.0358(4), uji = 0.845(10), uj2 = 2uji, and = 2 is illustrated in Fig. |l|. 
As we expected, the effective amplitudes of our four-parameter fit (solid symbols) 
tend to diverge with increasing of L, which shows that this set of critical exponents 
is false. One could object that, probably, the instability of the effective amplitudes 
is due to small errors in MC data. However, the amplitudes bi{L) and 62 (-L) of the 
more stable three-parameter fit (/ = 1, 2 in (|47|)) behave in a similar way (see empty 
symbols in Fig. ^). Moreover, the amplitude bi{L) of the two-parameter fit, shown 
by crosses, increases almost linearly at large enough L instead of the expected (in a 
case of correct exponents) saturation like bi{L) ~ 61 + const ■ . As regards the 
convergence in Fig. |l^ of the effective amplitudes at L ^ 00, it is possible only if 
both conditions are fulfilled, i. e., the exponents are correct and the MC data are 



accurate enough to ensure stable results. Thus, in any case, the analysis in Fig. |12 
provides rather convincing evidence that our exponents are the true ones, which 
by itself rules out the possibility that those proposed in |^5| could be correct. The 



results in Figs. 12 and 13 are affected insignificantly by a small inaccuracy of about 
10"'' in the estimated (3c value. 
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Figure 12: The effective amplitudes in Eq. (|^ 100 (Ina(L)-l) (triangles), 55i(L) (circles), 
b2{L) (squares), and 63 (L) (rhombs) evaluated by fitting the susceptibility data of 3D Ising 
model at criticality within the interval of sizes [L; 8L] with the critical exponents rj = 1/S 
and LJi — 1/2 of the GFD theory. 




Figure 13: The effective amplitudes in Eq. ( ^7|) evaluated by fitting the susceptibility data 
of 3D Ising model at criticality within the interval of sizes [L; 8L] with the critical exponents 
77 = 0.0358, uji = 0.845, uj2 = 2u;i, and lu^ = 2 proposed in ||2^. Sohd symbols show the 
four-parameter fit: 505i(L) (circles), 62 (-/j) (squares), and 63 (L) (rhombs); empty symbols 
show the three-parameter fit: 1005i(L) (circles) and 2762 (i) (squares); crosses represent 
the amplitude of the two-parameter fit, i. e., quantity 190 {bi{L) + 0.34). 
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11 Remarks about other numerical results 



There exists a large number of numerical results in the published literature not dis- 
cussed here and in |5|. A detailed review of these results is given in [31|. The cited 
there papers report results which disagree with the values of the critical exponents 
we have proposed. However, as regards the pure Monte Carlo study, we are quite 
confident that, just like in the actually discussed case of 3D Ising model, the in- 
crease of system syzes and/or use of higher-level approximations will lead to the 
conclusion that fits with our exponents are better than those with the conventional 
(RG) exponents. Particularly, a careful analysis of the effective exponents made 
in Sees. ^, ^, and |8| already has shown that the effective exponents deviate from 
the values predicted by the perturbative RG theory and converge more or less to 
those of the GFD theory at L ^ cxd. Together with the analysis of the experiment 
with superfluid '^He |p, we have presented totally 5 independent evidences of such 
a behavior. Besides, an important prediction of exceptionally our theory regard- 
ing the corrections to scaling is confirmed by the exact-algorithm transfer matrix 



calculations in Sec. 4.3 



There exists some background for the conventional claims in the published liter- 
ature that all the usual methods give consistent results which appear to be in a good 
agreement with the predictions of the perturbative RG theory. The perturbation 
expansions of the RG theory, as well as the techniques of high- and low-temperature 
series expansion are merely not rigorous extrapolation schemes which work not too 
close to criticality. As a result, these methods produce some pseudo or effective 
critical exponents which, however, often provide a good approximation just for the 
range of temperatures not too close to Tc (critical temperature) where these meth- 
ods make sense and, therefore, agree with each other. According to the finite-size 
scaling theory, tL^/'^ is a relevant scaling argument, so that not too small values of 
the reduced temperature t are related to not too large sizes L ~ t^^ . Therefore, one 
can understand that the MC results for finite systems often can be well matched 
to the conventional critical exponents proposed by high temperature (HT) and RG 
expansions. If, however, the level of MC analysis (i. e., the level of approximations 
used) is increased, then it turns out that the "conventional" critical exponents are 
not valid anymore, as it has been demonstrated in the current paper. It is because 
the "conventional" exponents are not the asymptotic exponents. Correct values of 
the asymptotic exponents have been found in considering suitable theoretical lim- 
its instead of formal expansions in terms of In k (at criticality, where k is the wave 
vector magnitude) or Int (approaching criticality) which are meaningless at A; — > 
and i — > 0. These formal expansions lie in the basis of the RG expansions for the 
critical exponents. One argue that Infc diverges weakly, therefore the expansions in 
powers of In k can be treated. This is a nonsense: any term like k~^ with A < can 
be formally expanded in terms of In A;, but therefore it does not become less diver- 
gent. Moreover, not only the powers of k but almost any function can be expanded 
in terms of In/c, therefore it is impossible to dechiper what is hiden behind such 
formal expansions in reality. This is a serious problem, since even an exponentially 
small correction (at A; ^ 0) can give a nonvanishing contribution to such a formal 
expansion (see examples in Sec. 2 of [^). The problem is not only formal: it has 
been proven in [S| that the assumption that the e-expansion works and provides 
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correct results at k ^ leads to an obvious contradiction in mathematics (cf. Sec. 2 
in (P). This fact alone cannot be compensated even by an infinite number of nu- 
merical evidences supporting the "conventional" critical exponents coming from the 
RG expansions. 

Our arguments, based on the current numerical analysis, are the following. First, 
the calculations by exact algorithms in Sec. 4.3 confirm our theoretical prediction, 
but not that of the perturbative RG theory. Second, we have proposed here a very 
sensitive method (i. e., a study of effective amplitudes) which allows to test the 
consistency of a given set of critical exponents with the MC data including several 
(in our case up to 3) corrections to scaling. We have applied this method to one 
of the recent and most accurate numerical data for the susceptibility in 3D Ising 
model, and have got a confirmation that our critical exponents are true. It would 
be not correct to doubt our results based on less sensitive methods and lower-level 
approximations . 

We prefer to rely just on the data of pure MC simulations becose of the following 
reasons. The so called Monte Carlo RG (MCRG) method is not free of assumptions 
related to approximate renormalization. We would like only to mention that the 
MCRG study in 1 33 ] of 3D Ising systems of the largest (to our knowledge) available 
in literature sizes, i. e. up to L = 256, has not revealed an excellent agreement with 
the usual predictions of the perturbative RG. In particular, an estimate uj ~ 0.7 has 
been obtained [^3| which is smaller than the usual (perturbative) RG value ~ 0.8, 
but still is larger than the exact value 0.5 predicted by the GFD theory. The high- 
temperature series cannot give more precise results than those extracted from the 
recent most accurate MC data, including the actual data of |27|, since these series 
diverge approaching the critical point. One approximates the divergent series by a 
ratio of two divergent series (Fade approximation) , but it is never proven that such a 
method converges to the exact result. Besides, the comparison to our calculations in 
2D Ising model via exact algorithms (Sec. ^ shows that the HT series analysis leads 
to misleading conclusions regarding such fine effects as corrections to scaling. These 
effects are relevant for 3D models. It is interesting to compare the MC and HT 
estimates of the critical point for the standard 3D Ising model, i. e., Pc — 0.2216545 
(MC) ^ and /3c = 0.221659 + 0.000002/ - 0.000005 (HT) [||]. It is clear that the 
MC value is more accurate: if we look in ||2^, where the estimation procedure is 
well illustrated, we can see that Pc is definitely smaller than 0.221659, and the error 
seems to be much smaller than the difference between both estimates 0.0000045. 
As we have mentioned already, our independent tests suggest that the error of the 
actual MC value is about 10"''. 



12 Conclusions 

Summarizing the present work we conclude the following: 

1. Critical exponents and corrections to scaling for different physical quantities 
have been discussed in framework of our Q recently developed GFD (grouping 
of Feynman diagrams) theory (Sec. |2|). 

2. Calculation of the two-point correlation function of 2D Ising model at the 
critical point has been made numerically by exact transfer matrix algorithms 
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(Sees. ^ and ^). The results for finite lattices including up to 800 spins have 
shown the existence of a nontrivial correction to scaling with a very small 
amplitude and exponent about 1/4 in agreement with the prediction of our 
GFD theory. No correction with the conventionally predicted exponent 4/3 
has been detected. 

3. The recently published Monte Carlo data for several three-dimensional lattice 
models have been reanalyzed. This analysis in Sees. || to |8| has shown that the 
effective critical exponents deviate from the values predicted by the pertur- 
bative RG theory and converge towards those of the GFD theory at L ^ oo. 
The same behavior has been observed in the experiment with superfluid ^He 
discussed in |^]. Totally, these are five independent evidences of such a behav- 
ior, suggesting that the above examples are not occasional or exceptional, but 
reflect a general rule. 

4. Different sets of critical exponents (one provided by GFD theory, another pro- 



posed in [25 1) predicted for the 3D Ising model have been tested by analyzing 



the effective amplitudes (Sec. ^ and While the usual fits of the suscep 



tibility data do not allow to show convincingly which of the discussed here 
sets of the critical exponents is better, this method strongly suggests that the 
conventional critical exponents rj = 0.0358(4) and uj = 0.845(10) ||25[ are false, 
whereas our (GFD) values ry = 1/8 and lo = 1/2 are true. 
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